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Proxemy Research is under contract to NASA to perform science research of volcanic plumes on Io 
and Mars. The following report constitutes delivery of Milestone Event #7 under NASA contract 
NASW-00013. 

TITLE: Io/Mars Quarterly Progress Report #101-04 

AUTHOR: Dr Lori S. Glaze 

1. Introduction 

Proxemy Research is under contract to NASA to perform science research of volcamc plumes 
on Io and Mars. Funding for this project began in May 2000 under contract NASW-00013 The 
project covered by this contract is comprised of two distinct tasks that were outlined in the original 
proposal entitled “Volcanic plumes on Io and Mars,” dated May 1999. The objective of the first task 
was to develop a model that constrains stochastic-ballistic effects of variable ejection velocities on 
areal concentrations of volcanic deposits on Io. The objective of the second task was to apply a 
buoyant plume rise model to explosive volcanic eruption columns released into an early Mars 
atmosphere 

2. Quarterly Report 

Tremendous progress was made in the last quarter of 2001 Over this three month period, 
efforts were directed toward three primary activities. (1) submission of a manuscript describing the 
limits of volcanic plume rise models to conditions on Mars [Attachment A], (2) preparation of a 
manuscript in collaboration with Ellen Stofan and others describing a statistical analysis of coronae 
on Venus, and (3) submission of a manuscript in collaboration with Stephen Baloga describing a 
probabilistic approach to pahoehoe lava flow emplacement [Attachment B], This last activity has not 
been funded directly by this contract, but is a collaborative PGG project. Thus it is reported here for 

completeness. 

Volcamc Plumes on Mars 

The second task of this contract is intended to focus on the application of plume rise models 
to explosive volcanic eruptions on Mars. The task considers both current and paleo-atmospheric 
conditions. Over the last two years, it has become clear that in the past terrestrial plume rise models 
have been applied to Mars without ensuring the validity of the models themselves. Our studies have 
shown that the convective rise models that have been used previously break down a few km above 
the eruption vent. As a result, other studies (e g., Wilson and Head, 1994; Hort and Weitz, 2001) 
have overestimated the heights to which volcanic eruption plumes can rise. 
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In early December a manuscript entitled “Volcanic plume heights on Mars: Limits of validity 
for convective models” was submitted to JGR/Planets. The manuscript, included here as Attachment 
A, describes where and why these convective rise models break down. We show that basic model 
assumptions are violated when (1) vertical velocities exceed the speed of sound, (2) radial expansion 
rates exceed the speed of sound, (3) radial expansion rates approach or exceed the vertical velocity, 
and (4) plume radius grossly exceeds the plume height. All of these assumptions are easily violated 
under current atmospheric conditions on Mars. The reason for the model breakdown is that the 
current Mars atmosphere is not of sufficient density to satisfy the conservation equations. It is likely 
that second order effects dominate the plume rise dynamics very rapidly. 

Coronae on Venus 

In late 2000, Proxemy Research was fortunate to employ Dr. Ellen Stofan. Throughout 2001, 
we have begun to develop several collaborations. One such collaboration has been to statistically 
analyze coronae on Venus that Dr. Stofan has previously classified in her corona database. Through 
this analysis, we are able to distinguish between various subpopulations of coronae. We expect to 
submit this manuscript within the next three months, and will attach the final version at that time. 

Pahoehoe lava flows 

I have been collaborating with Stephen Baloga on his PGG grant for many years. One of his 
primary tasks has been to explore a stochastic model for the emplacement of pahoehoe in the toe 
regime. This last quarter saw the submission of a manuscript to JGR/Solid Earth that we have been 
working on for quite some time. The manuscript, entitled “Pahoehoe transport as a correlated 
random walk”, is included here as Attachment B. The manuscript describes a correlated random walk 
model that can predict cross sectional profiles that are consistent with pahoehoe toe emplacement. 
Although this activity has not been firnded directly under this contract, it is included here for 
completeness. 
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Volcanic plume heights on Mars: 

Limits of validity for convective models 

Lori S. Glaze and Stephen M. Baloga, Proxemy Research, 20528 Farcroft Lane, 
Laytonsville, MD 20882, 301-313-0026, lori@proxemy.com 

Submitted to JGR/Planets, December 5, 2001 


Abstract. Previous studies have overestimated volcanic plume heights on Mars. In this work, we 
demonstrate that volcanic plume rise models, as currently formulated have only limited validity in any 
environment These limits are easily violated in the current Mars environment, and may also be 
violated for terrestrial and early Mars conditions. We indicate some of the shortcomings of the model 
with emphasis on the limited applicability to current Mars conditions. Specifically, basic model 
assumptions are violated when (1) vertical velocities exceed the speed of sound, (2) radial expansion 
rates exceed the speed of sound, (3) radial expansion rates approach or exceed the vertical velocity, 
(4) plume radius grossly exceeds plume height. All of these criteria are violated for the example given 
here. Solutions imply that the convective rise model is only valid to a height of ~10 km. The reason 
for the model breakdown is that the current Mars atmosphere is not of sufficient density to satisfy the 
conservation equations. It is likely that diffusion and other second order effects dominate the 
dynamics within the first few kilometers of rise. When the same criteria are applied to eruptions into 
a higher density early Mars atmosphere, we find that eruption rates higher than 1 .4 x 10 9 kg/s also 
violate model assumptions. This implies a maximum extent of ~65 km for convective plumes on early 
Mars. The estimated plume heights for both current and early Mars are significantly lower than those 
previously predicted in the literature. Therefore, global-scale distribution of ash seems implausible. 


Introduction 

Models of buoyant volcanic plumes 
have been applied to the Earth (e g.. Glaze et 
al. [1997]), Venus [Thornhill, 1993; Robinson 
et al., 1995; Glaze, 1999], Mars [Mouginis- 
Mark et al., 1988; Wilson and Head, 1994; 
Kusanagi and Matsui, 1998; Hort and Weitz, 
2001], and Triton [Soderblom et al., 1990; 
Kirk et al., 1990], These applications develop 
various inferences about the role of explosive 
volcanism in the transport and redistribution of 
volatiles, the character of eruption conditions, 
and the nature of surface deposits. These 
inferences are completely dependent on the 
validity of the underlying physics contained in 
the plume models. 


The role of explosive volcanism in the 
transport and redistribution of volatiles on 
Mars is particularly intriguing. The conjectures 
of a warm and wet paleoclimate [Pollack et al., 

1 987], possibly capable of sustaining rainfall 
[Craddock et al., 1998] and standing water 
[Parker and Currey, 2001], would favor 
buoyancy driven plumes. This, coupled with 
morphologic evidence that older volcanic 
centers were more explosive [Greeley and 
Spudis, 1981] may have meant a significant 
role for explosive volcanism in the injection of 
green-house gases into the atmosphere. 
Explosive volcanism under such conditions 
may have been important for the the 
redistribution of water on a regional basis as 
well as a factor in the global water budget. 
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Long runout pyroclastic flows, and other 
features on Mars may be the result of this type 
of paleovolcanism It has been proposed that 
layered terrain observed in high resolution 
MOC (Mars Orbiter Camera) images may be 
ash deposits from enormous volcanic plumes, 
that mantle broad regions [Malin and Edgett, 
2000], Furthermore, if explosive volcanism on 
Mars could produce plumes extending 
hundreds of kilometers into the atmosphere, 
these plumes could transport volatiles to the 
outer atmosphere, facilitating atmospheric 
release through Jeans escape [Chamberlain, 
1978], 

All of these important planetary science 
issues cause concern about the limits of 
validity for buoyant volcanic plume models 
and the progeny of inferences drawn from their 
applications. Many inferences, such as 
distribution of eruption products, are 
dependent on the maximum predicted plume 
height. Plume height estimates given in most 
terrestrial and planetary applications are based 
on conservation relationships established by 
Morton et al. [1956] for buoyant convective 
rise. The basis of the Morton formulation for 
convective rise is a system of first-order 
differential equations prescribing the 
conservation of mass, volume, and momentum. 
All aspects of turbulence are reflected in a 
single constant “entrainment” parameter. 
These models are intended to approximate 
solutions of the full Navier-Stokes fluid 
dynamics system. However, they contain 
many implicit assumptions and they have only 
been validated for laboratory, and small 
industrial plumes on Earth 

While there is broad agreement 
between theory and observations for volcanic 
plumes on Earth, it is important to remember 
that the simple first-order models are 
necessarily limited. This is due, in part, to the 
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fact that buoyant atmospheric plumes display 
physical processes that are extremely difficult 
to model for a wide range of circumstances 
and parameters. Volcanic plumes are 
turbulent, have horizontal boundaries that are 
free to respond to internal dynamics, interact 
intimately with the local ambient atmosphere, 
and contain a variety of gas, liquid and particle 
constituents. Also, from an empirical point of 
view, field measurements of the quantities 
needed to improve the theoretical models of 
turbulent atmospheric plumes are extremely 
difficult to obtain by comparison with 
laboratory or natural buoyancy (e g., industrial 
plumes) experiments with liquids. The 
application of these models to other more 
extreme environments, such as those found in 
planetary settings on Venus and Mars, 
exacerbates our need to understand exactly 
under what conditions the models are valid. 

Longstanding opinion in the literature 
is that explosive volcanic eruption plumes 
extend five times higher on Mars than analogs 
on Earth for similar eruption rates, and that 
column collapse occurs at heights twice those 
on Earth [Wilson and Head, 1994], Some 
estimates [Wilson and Head, 1994; Kusanagi 
and Matsui, 1998; Hort and Weitz, 2001] 
claim that explosive eruption plumes can climb 
to heights in excess of 100 km on Mars. This 
is more than twice as high as the Mount 
Pinatubo eruption plume, in an atmosphere 
that is considerably less dense. There are 
questions about the validity of convective rise 
models for volcanic plumes on Earth that have 
never been addressed. These questions 
become major concerns when models are 
applied to the Martian environment. 

It is our belief that, not only have basic 
assumptions of the models been violated in 
many applications found in the literature, but 
the buoyant plume models have been 
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extrapolated well-beyond their limits of 
validity. There are at least 7 areas of concern 
that limit the validity of plume rise models 
formulated in this way These areas of 
concern are as follows: 

1 ) Speed of sound violations: The models are 
totally inapplicable, despite what has appeared 
in the literature, when the speed of sound is 
exceeded (in the vertical or radial directions), 
or the atmosphere is too thin to permit 
definition of the speed of sound [Stothers, 
1989], 

2) Radial expansion: The intent of these 
models is to describe a plume that is moving 
primarily upward, not outward. As a result, 
the rate of radial expansion should be small 
compared to the vertical velocity and the 
overall plume radius should be small compared 
to the plume height. 

3) Ambient air availability: Turbulent 

convection models implicitly assume that there 
is enough ambient air for the plume to behave 
as a continuum. The models also assume that 
the corresponding dynamic impact of the 
plume on the atmosphere is negligible 

4) Instant mixing: Most volcanic applications 
of the Morton approach have considered only 
the bulk properties of the plume. This carries 
with it the implicit assumption that all changes 
in pressure, temperature and density within a 
control volume are felt instantly and 
simultaneously throughout each horizontal 
control volume 

5) Confluent particle and gas flow: Particles 
provide the majority of the heat that drives 
buoyancy. When we assume that particles and 
gas flow confluently at all times, it is implied 
that there must be enough gas pressure to 
suspend the particles and only enough 


turbulence to do this without ejecting them. In 
addition, the heat contained in the particles 
must be made available to the gas. The time 
required for conductive heating of the gas 
must be negligible. 

6) Pressure balance: One of the most basic 
simplifying assumptions is that the pressure 
within the plume is the same as without. 
Exceeding the speed of sound is an obvious 
violation of the pressure balance assumption, 
but there are other indicators that the pressure 
balance assumption is violated even in lower 
velocity cases, including rotating columns 
observed in firestorms and industrial plumes. 

7) Diffusion and other second order effects: 
Many second-order effects accounted for in 
the complete Navier-Stokes treatment (such as 
diffusion of mass, turbulent eddies, momentum 
and heat) are ignored in the Morton system. In 
this reduction, all the second-order effects 
have been compressed into the entrainment 
constant, a. It is not clear when the ignored 
second order effects become important, or if 
the entrainment parameter is indeed constant 
in different turbulent environments (eg.. 
Mars). 

Thus, the question we hope to answer 
with this paper is, “What are the limits of 
validity for the Morton-type convective rise 
formulation for conditions on Mars?”. We 
begin with a brief review of the approach that 
has been the most common method for 
describing volcanic plumes, including the basic 
assumptions in the original Morton 
simplification. We then analyze to what 
extent these assumptions are valid by 
addressing the concerns listed above. We 
focus on current atmospheric conditions found 
on Mars, as well as the Earth and early Mars. 
Finally, we attempt to answer the questions: 
how high could a volcanic plume realistically 
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rise in the current Mars atmosphere, and how ^ 

high would the same plume go in a ‘paleo’ j = ^ a P aB ur ( 1 ) 

atmosphere? 


Convective Rise Models 

For completeness, we provide here 
some background and basic information 
concerning the convective rise models that are 
referred to throughout this paper. Over the 
past two decades, several models describing 
the convective rise of terrestrial volcanic 
plumes using a system of first-order 
conservation equations have appeared in the 
literature [Wilson et al., 1978; Sparks and 
Wilson, 1982; Sparks, 1986; Woods, 1988, 
1993; Glaze and Baloga, 1996; Glaze et al , 
1997], These models are based on the original 
work of Morton et al [1956], Morton et al. 
attempted to approximate the more complex 
Navier Stokes approach to modeling small- 
scale buoyant flows. It is interesting to note, 
however, that as yet no one has ever shown a 
direct link between the time-honored Navier- 
Stokes and the Morton system. Morton et al. 
approximated the turbulent entrainment of 
ambient fluid into a convectively rising plume 
by using an empirical parameter, referred to 
here as a. This was originally a very 
innovative approach that side-stepped the 
complex modeling of both a turbulent plume 
and the motion of the ambient atmosphere 
with a free-boundary between them. 

The basic set of conservation equations 
to be examined in this study are taken from 
Glaze et al. [1997], hereafter referred to as the 
model. Other models may contain slight 
differences, but the basic assumptions of 
Morton et al. [1956] dominate all results. The 
basis of this approach is the expression for the 
conservation of mass, given by 


Please see the Notation table for definitions of 
all variables. The entrainment constant, a, 
appears explicitly in (1), indicating that the 
mass of ambient fluid entrained is directly 
proportional to the upward rise velocity of the 
plume itself. 

The momentum conservation equation 
corresponding to (1) that prescribes the 
buoyancy relationship with the ambient fluid 
can be written as 

■^{pB u2r2 ] = g{PaB-pB) r2 ( 2 ) 

The thermal energy conservation is 
) B «r 2 C fi 0] = 2ap ^urC^T- PaB u r 2 g$ d Q) 

where the first and second terms on the right 
hand side describe how the bulk plume 
temperature changes due to the entrainment of 
ambient fluid at a different temperature and 
adiabatic expansion, respectively. 

This simplified approach to modeling 
buoyantly convecting plumes was stimulated 
by a desire to understand the dynamics of 
small industrial plumes. To this end, detailed 
validation studies were conducted (e g., Briggs 
[1969]) to ensure that the models could 
adequately describe the behavior of convecting 
plumes in the lower troposphere. Similar 
models have also been applied to the dynamics 
of firestorms that are sources of enormous 
heat [Carrier and Fendell, 1983], Wilson et al. 
[1978] and Settle [1978] compared model 
results to observations of volcanic plume rise 
and concluded that the Morton-type model of 
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convective rise was broadly applicable to 
plinian style volcanism on Earth. Sparks et al. 
[1997] also present data that support the use 
of convective plume models, such as that 
shown above, to describe volcanic eruption 
columns on Earth. These volcanic comparisons 
focused on smaller, tropospheric plumes. 
Thus there are still unresolved issues regarding 
the broad applicability of the Morton style 
models to volcanic plumes that penetrate the 
stratosphere. Furthermore, these comparative 
analyses are not strictly validation studies, as 
the true eruptive parameters are always 
unknown. 

There are many fundamental 
assumptions inherent to all the first order 
models for convective rise [Stothers, 1989], 
The most basic of these include the 
assumptions that the (1) pressure within the 
plume is the same as ambient for all altitudes 
(i.e, plume is ‘pressure balanced’), (2) plume 
material is transported primarily upward (i.e., 
the vertical velocity is large relative to the 
radial expansion), (3) particles are able to 
transfer heat instantaneously to the plume gas 
(heat conduction within the particle is not 
important), (4) density difference between the 
plume and ambient fluid is “small” (although 
neither the absolute difference, nor the 
difference relative to other parameters, 
necessary for model validity have ever been 
quantified), (5) particles and gas move 
together (e g., no drag), and (6) ambient 
atmosphere is of sufficient density that the 
plume can be treated as a continuum An 
implicit corollary to point (1) is that the 
vertical velocity of the bulk plume must be less 
than the local speed of sound. 

In applying such models to volcanic 
eruption plumes on Mars, point (6) becomes 
of critical importance The very core of the 
model described above is that the eruption 
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column must convect, not merely rise 
buoyantly. Buoyancy is simply the result of 
density differences in fluids. For example, a 
balloon can rise buoyantly, but it cannot 
convect. The process of convection requires 
the entrainment of ambient fluid and, thus, 
sufficient ambient fluid must be available for 
convection to occur at all. The existence of 
sufficient ambient fluid is a subtle requirement. 
But, without ambient fluid, there is nothing to 
“entrain”. 

Comparison of Atmospheric 
Conditions 

We now compare several attributes of 
the current Mars atmosphere with the Earth 
and plausible paleo-atmospheric conditions on 
Mars. We then present several examples 
illustrating how and when the basic 
assumptions above break down. 

The convective rise models described 
above are critically dependent on the basic 
assumption that there is sufficient atmosphere 
available for entrainment. We know that at 
some limit in the absence of sufficient 
atmosphere, the dynamics of the ejected 
material must be controlled by ballistic 
behavior, such as that observed on Io [Glaze 
and Baloga, 2000], Where and if the transition 
occurs is very difficult to ascertain, but the 
point at which the definition of the speed of 
sound breaks down is certainly an upper limit 
for a continuum model. Before we can begin 
looking at the model, however, we must first 
characterize the atmosphere in which we are 
trying to convect. In this section we define the 
basic atmospheric parameters that are 
important to the assumptions of the convective 
models and provide comparisons with our own 
atmosphere on Earth. 
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Current atmospheric conditions on 
Mars are significantly different than on Earth. 
The gravitational force holding the atmosphere 
is less, the composition is different (primarily 
C0 2 on Mars vs N 2 and 0 2 on Earth), the 
surface pressure is two orders of magnitude 
less, and the atmosphere is considerably 
cooler. While the pressure of the paleo-Mars 
atmosphere may have been closer to current 
conditions on Earth, the composition, 
temperature and gravity still strongly influence 
convective plume rise. 

For comparative purposes, we use the 
US Standard Atmosphere for Earth (eg., 
Handbook of Chemistry and Physics [1955]). 
For Mars, we use parameters derived by Seif 
and Kirk [1977], More recent estimates of 
atmospheric conditions on Mars derived from 
the Mars Pathfinder descent [Magalhaes et al., 
1999] are similar to the Viking-derived data. 
For paleo-conditions on Mars, we allow the 
surface pressure to be the same as Earth, 
keeping everything else the same as the current 
Mars atmosphere. 

The temperature profiles used in this 
comparison and in the modeling examples 
below are shown in Figure 1 . The temperature 
profile for Mars has been simplified for 
modeling purposes. Pressure profiles, shown 
in Figure 2, are derived from the expression 
P(z) = P(0) exp(-z/H), where P(0) is the 
surface pressure, and H is the atmospheric 
scale height. Table 2 summarizes these and 
other constants assumed for each of the model 
atmospheres. 

The density profile for each 
atmosphere is strongly controlled by the 
pressure. To find the density at each step, we 
have used the ideal gas relationship p^z) = 
P(z)/R aB T(z). To get a better feel for the 
density of each atmosphere, we have 


calculated the mean free path as a function of 
altitude. The mean free path is the average 
distance that a molecule travels before 
colliding with another molecule. Most 
introductory physics texts (e.g.. Sears and 
Zemansky [1967]) discuss the mean free path 
to some extent. The mean free path, L(z), can 
be estimated using the expression 


L(z) = 


sin(7i/4) 
o n(z) 


(4) 


where n(z) is the molecular density . 

R ob mw 7 ’( z ) 


(5) 


and N is Avogadro’s Number (N = 6.0225 x 
10 23 ). For consistency with the Handbook of 
Chemistry and Physics, we have used o = 
4.165 x 10' 15 cm 2 in (4). 

Figure 3 is a plot of the mean free path 
as a function of altitude on both Earth and 
Mars. For illustrative purposes, we have 
highlighted a few interesting things as well. 
First, we have indicated the maximum 
estimated altitude of the Pinatubo column from 
the 1991 explosive eruption of ~ 35 - 40 km 
[Self et al., 1994], While other eruptions 
produced larger volumes than Pinatubo, there 
is no evidence that plumes from these 
eruptions would have extended much further 
into the atmosphere. The maximum plume 
height during the Krakatau eruption of 1883 
has been estimated at ~45 km [Self, 1992], 
The eruption plume from Toba (73,500 bp), 
the largest known explosive volcanic event of 
the late Quaternary, has been estimated at only 
27 - 37 km [Rampino and Self, 1992], The 
seeming coincidence that all these plumes 
seem to have similar maximum heights 
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suggests that there may be some empirical 
terrestrial limit around 40 km that is not 
explained by current plume models 

The mean free path at this empirical 
terrestrial limit is very similar to the mean free 
path at the surface on Mars. The range of 
mean free paths where we see convecting 
plumes on Earth is between 10‘ 3 and 10' 5 cm. 
Out of interest, we have also highlighted the 
fact that the average distance between 
molecular collisions at altitudes where 
convective plume rise has been predicted on 
Mars is 1 - 100 meters! Convectively rising 
plumes must quickly attain an internal density 
that is similar in magnitude to the ambient 
atmosphere. This means that at these extreme 
altitudes, the volcanic plume must also have a 
mean free path equivalent to the atmosphere at 
that height. Intuitively, it seems unlikely that 
particles that must travel even a meter before 
a collision are capable of entraining anything. 

We have also indicated that even the 
U2 aircraft can only fly as high as about 22 
km. Beyond this altitude, there is not enough 
atmosphere to create aerodynamic lift. 
Eyewitness accounts of atomic bomb tests 
indicate that these powerful blasts cannot 
produce plumes that exceed 9 - 12 km. 
Furthermore, the firestorm of Hamburg, 
Germany, in 1943 attained temperatures in 
excess of 1000 C for a prolonged period, 
perhaps as much as 6 hrs Some estimates 
suggest the temperature may have reached 
3000 C, yet the plume only attained a height of 
9-12 km kilometers [Carrier and Fendell, 
1983] 

The most basic assumptions of the 
convective rise models is that of balanced 
pressure inside and outside the plume. One of 
the most obvious ways to violate this 
assumption is to have vertical rise or radial 
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expansion velocities that exceed the speed of 
sound. Haliday and Resnick [1977] give the 
speed of sound as 



where y is the ratio of specific heats (Table 2). 
For N 2 and O a which make up 98.9% of 
Earth’s atmosphere, y = 1.4. For a primarily 
C0 2 atmosphere, y = 1.3. Figure 4 compares 
the speed of sound on Earth and Mars. 
According to the Handbook of Chemistry and 
Physics, the concept of the speed of sound in 
Earth’s atmosphere loses its applicability at 
about 90 km where the mean free path of air 
molecules approaches the wavelengths of 
sound waves (~ 4 cm). This loss of definition 
occurs at about 70 km on Mars. 

The balanced pressure assumption can 
also be violated even when plume velocities 
are less than the speed of sound. Studies of 
firestorms indicate that the basic pressure 
balance assumption breaks down with the 
onset of a rotational motion in the column that 
is not accounted for by Morton types of 
formulations. This rotation causes a pressure 
gradient and may reduce the ability to entrain 
ambient air (i.e., a) by as much as two orders 
of magnitude. The validity of Morton 
formulations has been tested only up to the 
Long Beach firestorm, which attained a 
maximum plume height of only several 
kilometers [Carrier and Fendell, 1983], 

Validity Issues 

The discussion above regarding 
ambient atmospheric conditions suggests that 
there may be some point in Mars’ atmosphere 
where the ambient density is not sufficient for 
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plumes to rise convectively. Certainly, there is 
cause for further investigation. The crux of 
the problem is where to establish the limit of 
validity for the system of ordinary differential 
equations in (1) - (3). Here, we attempt to 
establish where the model breaks down. For 
the purpose of the following discussions, we 
will use the basic boundary conditions in Table 
3. These are typical boundary conditions that 
have been used for explosive eruptions on 
Mars that can be found in the literature (e g., 
Hort and Weitz [2001]). For comparison, we 
have shown the parameters for an eruption on 
Earth, and into a more dense early Mars 
atmosphere, with equivalent mass eruption 
rates. Unless otherwise indicated, the volatile 
is assumed to be water vapor. 

It is critical to this discussion to 
understand that the model cannot predict the 
maximum plume height directly because the 
system of equations (1) - (3) have a singularity 
at u = 0. As a result of this singularity, both 
the radius and velocity terms change very 
rapidly near the top [Glaze and Baloga, 1996], 
Previous terrestrial studies have estimated the 
maximum predicted plume height on Earth by 
taking the asymptotic solution near this 
boundary layer. In our previous studies [Glaze 
and Baloga, 1996; Glaze et al., 1997], this 
asymptotic solution has been taken as the point 
at which the upward velocity falls below some 
set value (usually 10 m/s for large eruptions). 
On Earth, this cutoff seems to give a 
reasonable estimate of the limit of validity for 
the convective rise model. The predicted 
heights shown in Table 3 correspond to the 
points at which the velocity falls below 10 m/s 
as a default case On Mars, however, this 
criterion seems insufficient. As can be seen 
from Table 3, the upward velocity does not 
drop below 10 m/s until 133 km above the 
vent. This maximum height is consistent with 
other predicted plume heights cited in the 


literature for similar mass eruption rates 
[Wilson and Head, 1994; Kusanagi and 
Matsui, 1998; Hort and Weitz, 2001], The 
following discussions indicate a variety of 
other considerations that make applications of 
the model invalid well below this altitude. 

Io Comparison: 

We know it is absurd to assume 
convection on Io where we have visual 
evidence of ballistic dynamics (e g., Strom et 
al. [1981] for early Voyager reports; and 
Glaze and Baloga [2000] for modeling studies 
of ballistic dynamics). However, the 
convective model discussed here cannot make 
this distinction on its own. If we provide it 
with boundary conditions, it will attempt to 
predict convective plume rise, even under 
absurd conditions. 

As an illustration of this, we have run 
the model for an average Io atmosphere 
(understanding that the atmosphere on Io 
changes dramatically between day and night). 
For modeling purposes, we have assumed an 
average surface temperature of 110 K that 
warms slightly with altitude, a surface pressure 
of 1.7 x 10-* atm (.0017 Pa), and an 
atmospheric composition of S0 2 . These 
constraints result in an atmospheric density 10 
million times less than Earth. For boundary 
conditions similar to those on Mars in Table 3 
(except gravitational force of 1.8 m s' 2 ), and a 
magmatic volatile composition that is 50% 
S0 2 and 50% H z O, the model blindly predicts 
that a plume can convect to heights in excess 
of 1 50 km on Io. 

While we know that explosive 
eruptions on Io do produce plumes up to 300 
km high [McEwen and Soderblom, 1983], the 
conditions described above would only 
generate a ~15 km ballistic plume. These 
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plumes are clearly driven by very different 
mechanisms than those assumed by the model. 
For the Mars case then, the primary question 
is where do the basic assumptions of the 
convective rise model break down? 

Speed of Sound: 

One of the most basic assumptions of 
the convective rise model is that the pressure 
within the plume is always identically equal to 
the ambient pressure at that altitude. If the 
plume is moving slowly enough and if the 
radius is small enough, this is a plausible 
assumption. In other words, as long as the 
time it takes for the effects of the ambient 
pressure to propagate throughout the entirety 
of the plume volume is small compared to each 
time step, this assumption should be 
reasonable. While it is difficult to identify a 
specific cutoff where this assumption is 
violated, it is clearly violated at any point when 
the vertical rise velocity exceeds the ambient 
speed of sound 

It has been postulated that initial 
eruption velocities on Mars may well exceed 
the speed of sound [Wilson and Head, 1 994; 
Hort and Weitz, 2001], However, convective 
rise models should never use such velocities as 
boundary conditions. As a note, the 
convective models are not particularly 
sensitive to the eruption velocity. For 
example, increasing the initial velocity for 
Mars in Table 3 to u Q = 500 m/s, the 10 m/s 
plume cutoff predicts a maximum height of 
141 km (a 6% difference). For this study, we 
only examine results for vent velocities that are 
below the estimated speed of sound at the 
surface (-236 m/s). 

Figure 5 shows the velocity profile 
predicted by the convective model using the 
current Mars boundary conditions defined in 
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Table 3. The general shape of the velocity 
profile is similar to that predicted for large 
volcanic eruptions on Earth. For large 
eruptions, the bulk plume density just above 
the vent can be significantly greater than the 
ambient density. Thus, the plume is not 
buoyant and the column material will slow 
down as it rises. A convecting plume, 
however, is simultaneously entraining ambient 
atmosphere which expands as it is heated by 
the plume. If the eruption has enough energy, 
the plume will become buoyant before it loses 
all of its initial momentum. 

The unusual feature of the modeled 
velocity profile for Mars is that the 
convectively driven velocity far exceeds the 
initial eruption velocity. In addition, the 
upward velocity exceeds even the ambient 
speed of sound at an altitude of ~34 km. 
Thus, the convective rise model cannot be 
valid above this point for two reasons. First, 
the model in no way accounts for the 
atmospheric shocks caused by such velocities. 
Second, continuing past this point certainly 
violates the basic pressure balance assumption 
upon which the convective rise model is based. 

Radial Expansion: 

Ignoring the vertical velocity 
arguments made above, we need to recall that 
all the variables in the system of equations (1) 
- (3) are solved simultaneously. This means 
that if one wants to accept the predicted plume 
height given in Table 3, one must be willing to 
accept ALL the solutions for all the variables 
given in (1), (2) and (3). In particular, there 
are several anomalous features of the predicted 
plume radius that must be accounted for in 
some way. 

The first, and possibly most alarming 
anomaly, is the predicted radius itself. Figure 
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6 shows the radial solution that corresponds 
identically to the velocity solution in Figure 5. 
If one accepts the final plume height at the 
point where u < 10 m/s, then one must also 
accept a plume radius at that point that is in 
excess of 44,000 km! The circumference of 
Mars at an altitude of 133 km is ~ 22,000 km, 
so this would imply that the entire plume top 
could wrap around Mars four times. 

The dramatic increase in the radius 
implies entrainment of an enormous volume of 
ambient atmosphere. Aside from the sheer 
volume of entrained air, the model completely 
ignores the simultaneous dynamics that would 
be required for an atmosphere to 
accommodate this scale of ingestion. 

We can immediately see from visual 
inspection, at the scale of Figure 6, that the 
radial variable is clearly in a boundary layer 
regime above 110 - 120 km where the 
convective model cannot be valid. But where 
does the radial solution go wrong? Figure 7, 
shows the radius solution in the altitude range 
up to about 1 10 km. However, the radius is 
already changing very rapidly by about 80 km. 
This implies that the model may break down 
even before this point. 

To understand what is happening with 
the radius solution, we have estimated the 
radial expansion velocity (dr/dt) as a function 
of altitude. To do this, we simply difference 
the predicted radius on either side of each step 
in z. We then divide this Ar by the time it 
takes the control volume to travel the 
integration step size (5 m was used in all 
calculations) using the predicted velocity at 
that step. 

Figure 8 shows the radial expansion 
velocity along with the vertical velocity and 
speed of sound for comparison. We see 
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immediately, that the radial expansion velocity 
exceeds the ambient speed of sound even 
before the vertical velocity, at a height of ~27 
km. It can also be seen from Figure 8 that 
there is a complex relationship between the 
vertical and radial velocities. One of the 
implicit assumptions of the convective rise 
model is that the primary movement of the 
plume is upward, not outward (i.e., the radial 
velocity should be small compared to the 
vertical). However, the radial velocity far 
exceeds the vertical velocity at altitudes above 
~22 km. In addition, there is some question as 
to the validity of the initial assumption at 
heights above ~9 km where the two 
components almost cross. 

The driver that is causing the radial 
solution to break down is not immediately 
obvious. Figure 9, however, does shed some 
light on the issue. In this figure, both the bulk 
plume density and the ambient atmospheric 
densities are plotted as functions of altitude. It 
can be seen from Figure 9 that the plume 
becomes buoyant (less dense than ambient) at 
about 10 km. Above that point, the density 
difference is extremely small (10‘ 4 kg/m 3 at 10 
km to 10' 5 kg/m 3 at 50 km). From the form of 
(2) it can be seen that in order to conserve 
momentum when (p^ - p B ) is very small, r 
must become quite large to compensate. The 
fact that this density difference is so small 
raises several important questions. The first is 
what happens when there are small 
fluctuations or deviations from our 
assumptions in the atmospheric conditions? 
Glaze and Baloga [1996] showed that 
systematic variation of 2%/km in the cooling 
rate of a converting plume can reduce the 
predicted plume height by 1/3. Furthermore, 
the model assumes that any changes in plume 
temperature and density are felt throughout 
the plume volume instantaneously. For plume 
radii observed on Earth that exceed several 
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kilometers, this requires some stretch of the 
imagination. For plumes with predicted radii 
exceeding 100's of km on Mars, it seems 
ridiculous. 

Discussion 

Based on the discussions above, we 
believe that volcanic plumes cannot rise 
convectively to excessive heights in the current 
Mars atmosphere. For the example given, 
convection quickly breaks down within the 
first 10 km. This conclusion is in stark 
contrast to other studies that predict 
convective plumes in excess of 100 km. While 
volcanic plumes may indeed go higher than 1 0 
km, there must be some other mechanism that 
dominates the transport dynamics. 

One possible reason for the breakdown 
of the convection model is most likely in the 
exclusion of 2nd order effects, such as 
diffusion of momentum, turbulence, or energy. 
It may be that for conditions on Mars, these 
effects play a much more significant role. 
Diffusion and rotation serve as an energy sink 
and may cause the bulk upward velocity of the 
plume to decrease at a rate greater than 
predicted by the simplified approach. With the 
plumes on Mars we have found a unique case 
where the density difference between the 
plume and ambient is extremely small, and yet 
the plume density can never quite “catch up” 
to the atmosphere. Another possibility is the 
inability to maintain the confluence of gas and 
particle motions. As the particles carry the 
bulk of the heat content, this may dramatically 
change the internal dynamics of the plume. 

Unfortunately, none of this discussion 
helps us to estimate how high material erupted 
explosively on Mars may be transported 
realistically. We can, however, constrain the 


minimum plausible limits. The first constraint 
that we can place on plume rise is a minimum 
estimate based on a purely ballistic vertical 
trajectory. The maximum height attained by a 
ballistic mass is solely dependent on the 
velocity at the vent and gravity, H = u e 2 /2g. 
For a vent velocity of 230 m/s, we can safely 
say that a plume would have to be a minimum 
of 7 km high. 

Now, one possibility is to apply the 
convective model up to about 9 km (where we 
believe we are not violating any basic 
assumptions), and then allow the plume to 
continue upward ballistically. At 9 km, the 
plume density is still greater than the ambient 
(i.e., not yet buoyant) but, the plume would 
continue to rise ballistically as a result of its 
momentum. The model predicts that the plume 
has a bulk vertical velocity of ~106 m/s at 9 
km. Thus, if no more entrainment/convection 
occurs, the plume would rise an additional 1 5 
km, for a total of 10.5 km. 

One interesting alternative to consider 
is the possibility that there is never enough 
atmosphere for the plume to convect, but that 
the plume material could still rise due to 
momentum and buoyancy. This is equivalent 
to setting the entrainment constant equal to 
0.0 (no entrainment). For the initial conditions 
given in Table 3 and a = 0, a plume would rise 
to a height of about 8 km. 

All of these activities lead to the 
inference that explosive eruptions on Mars 
could plausibly inject material to altitudes 
between 10, or perhaps 20 km, under current 
conditions. The impact of this conclusion is 
that these eruptions would not be capable of 
distributing material over significant distances. 

ft is also of interest to consider 
eruptions into the early Mars atmosphere. 
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Assuming a surface pressure of the same order 
of magnitude as that currently found on Earth, 
the density of the early Mars atmosphere 
would have been more dense than our own 
The higher density results from bulk gas 
constants and atmospheric temperatures that 
are less than Earth. The higher density 
atmosphere and lower gravity should both 
result in conditions conducive to convective 
rise 

But, just how high would we expect 
plumes to go, and would any of the concerns 
discussed above become issues? In Table 3 we 
show the predicted maximum height for a 
plume with the same mass eruption rate as the 
example in the current Mars atmosphere. The 
convective model predicts that this plume 
would rise to a height of 3 1 km, about 50% 
higher than the same plume on Earth. This 
result is consistent with intuition that plumes 
should go somewhat higher in a lower gravity, 
higher density environment. The initial 
boundary conditions for this example are also 
within reason, with a vent radius of 1 00 m, and 
eruption velocities below the speed of sound. 
Close inspection of the velocity and radial 
solution indicates that the violations noted 
above do not occur under these conditions. 
Thus, it seems reasonable to model this event 
as a convecting plume. 

We have also explored the possibility 
of greater eruption rates to determine if there 
is a limit to how high plumes could go in an 
early Mars atmosphere. The convective rise 
model does, indeed, have a limit of 1.4 x 10 9 
kg/s. Mass eruption rates exceeding this limit 
produce plumes that violate all the criteria 
discussed for eruptions into the current 
atmosphere. This mass eruption rate is only 1 
order of magnitude greater than the examples 
give in Table 3, and is equivalent to a vent size 
of 300 m (with everything else held the same). 
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The maximum predicted plume height resulting 
from this mass eruption rate is approximately 
65 km. Again, this apparent limit is about 
50% higher than the empirical limit observed 
on Earth. 

The limit of 65 km under early Mars 
conditions, while somewhat higher than 
terrestrial plumes, is probably not sufficient to 
distribute volcanic ash globally. While 
aerosols are known to have stayed aloft for 
months to years following the Pinatubo 
[Trepte, 1993] and El Chichon eruptions 
[Pollack et al., 1983], ash particles are 
believed to have fallen out over a matter of 
hours to days. The most important implication 
of this conclusion is that it seems unlikely that 
volcanoes could be responsible for broad scale 
mantling and the layered terrain observed in 
MOC images. Volcanoes could have, 
however, been major suppliers of volatiles and 
aerosols into the atmosphere as well as a 
mechanism for volatile distribution processes. 

Conclusions 

Based on our analyses, maximum 
plume heights on Mars extend only a few tens 
of kilometers and are comparable to those on 
Earth. Volcanic plumes that rise to heights in 
excess of 100 km in Mars’ current atmosphere 
are highly improbable. The low density of the 
current atmosphere simply cannot support 
convection. The extreme maximum plume 
heights predicted by other investigators have 
resulted from a breakdown of the physics 
model used to describe the plume dynamics. 
We have demonstrated the breakdown of 
model validity in a variety of ways: 

1 ) Speed of sound violations: The models are 
totally inapplicable when the speed of sound is 
exceeded or the atmosphere is too thin to 
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permit definition of the speed of sound. For 
the example of a moderately large eruption 
discussed, we have shown that both the 
vertical and radial velocity components violate 
this assumption at 34 km and 27 km, 
respectively. Thus the model is most certainly 
invalid at heights above 27 km where no 
account is made for atmospheric shocks or the 
fact that the plume is no longer pressure 
balanced. The assumption of pressure balance 
may even be violated long before the velocity 
components exceed the speed of sound (on 
Earth as well as Mars). 

2) Radial expansion: One of the most basic, 
implicit assumptions of convective rise models 
is that plumes rise primarily upward, not 
outward. However, the radial expansion 
velocity of the example plume is of similar 
magnitude to the vertical velocity by about 9 
km, and actually exceeds the vertical velocity 
at a height of 22 km. In addition, were we to 
trust the model results above these heights, the 
radius of the plume quickly exceeds the plume 
height by orders of magnitude. At a minimum, 
this excessive expansion indicates that one 
should be focusing on the radial dynamics, as 
opposed to vertical. Moreover, the entrained 
ambient air must be mixed throughout each 
horizontal cross section of the plume, within 
the time it takes to travel each time step. As 
the plume continues to expand, this horizontal 
mixing assumption becomes increasingly less 
plausible. 

3) The availability of ambient air: The reason 

for excessive radial expansion is that there is 

not sufficient ambient air available to the 

plume for entrainment. The density of the 

atmosphere is so low, that the plume must 

expand large distances to find sufficient 

atmospheric mass to satisfy the conservation 
equation requirements. However, the model 
in no way accounts for the affect of the plume 


on the dynamics of the atmosphere. At a 
minimum, the atmosphere requires separate 
modeling. 

It may well be that volcanoes can inject 
materials to extreme heights, but it seems 
highly improbable that convection could 
possibly be capable. We have shown that for 
mass eruption rates equivalent to moderately 
large eruptions on Earth, an eruption plume on 
Mars could only convect as high as about 9 
km. A combination of convection and ballistic 
transport results in plume heights of only 10 
km. It is possible that material produced in 
explosive volcanic eruptions could go 
somewhat higher, but it would have to be 
through some other, mechanism that has not 
yet been modeled. 

An early Mars atmosphere was much 
more capable of sustaining convecting plumes. 
Consistent with results of other researchers, an 
eruption of similar magnitude to the current 
Mars example would climb to a height of 3 1 
km, ~50% higher than the equivalent plume on 
Earth. However, there appears to be a limit to 
the validity of the models on early Mars as 
well. This limit of validity occurs at an 
eruption rate of 1.4 x 10 3 * * * * * 9 kg/s (1 order of 
magnitude larger than the 3 1 km high plume) 
and produces plumes with a maximum height 
of 65 km. Eruption conditions exceeding this 
mass eruption rate violate model assumptions 
in the same ways as described above. 

The net result of these studies is that 
explosive volcanic eruption plumes can not 
rise as high as previously projected in the 
literature. Consequently, injection heights for 
ash and volatiles have been overestimated, as 
well as the fall heights of pyroclastic flows. 
Layered deposits and mantles are currently 
being recognized in the high resolution Mars 
images. Conjectures often arise about volcanic 
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airfall or pyroclastic flows from distant 
sources. Our results indicate that these 
sources should not be much farther than 
analogous deposits on Earth Whether the 
current atmosphere or paleo-atmo sphere is 
considered, injection of only a few tens of 
kilometers into the atmosphere need be 
considered. Thus, this localization of volcanic 
effluents implies a local to regional emphasis 
for volatile transport processes in contrast to 
large-scale, global atmospheric circulation. 

Our recommendation for improving the 
buoyant plume rise models is to move toward 
a more generalized form the Navier-Stokes 
equations, which are second-order, 
supplemented by appropriate constraining 
equations (atmosphere, equations of state, 
etc). It is most essential to begin addressing 
the diffusion of mass, turbulent eddies, 
momentum, and heat in the transitional regime, 
where differences between atmospheric 
properties (velocities, eddies, heat) are really 
second-order effects. For the unique 
conditions on Mars, it may well be that these 
second-order effects become important within 
the first 10 - 20 km of rise. 
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Table 1. Notation 


Variable 

Definition 

c B 

Bulk specific heat for plume (J K" 1 kg' 1 ) 

CaB 

Bulk atmospheric specific heat (J K' 1 kg' 1 ) 

g 

Acceleration due to gravity (m s' 2 ) 

H 

Maximum plume height (m) 

L(z) 

Mean free path (m) 

mw 

Molecular weight (kg/mol) 

N 

Avogadro’s number (6.0225 x 10 23 ) 

n(z) 

Molecular density (molecules m' 3 ) 

P(z) 

Pressure (Pa) 

r 

Plume radius (m) 

^aB 

Bulk atmospheric gas constant (J K' 1 kg' 1 ) 

T( z ) 

Atmospheric temperature (K) 

u 

Vertical plume velocity (m s' 1 ) 

V 

Speed of sound (m s" 1 ) 

z 

Vertical coordinate (m) 

a 

Entrainment parameter (= .09) 

Y 

Ratio of specific heats 

0 

Bulk plume temperature (K) 

Pb 

Atmospheric density (kg m' 3 ) 

PaB 

Bulk plume density (kg m' 3 ) 

0 

Mean molecular cross sectional area (m 2 ) 

4>d 

Fraction of control volume occupied by entrained gas 

4>v 

Fraction of control volume occupied by magmatic gas 
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Table 2. Summary of atmospheric parameters. 
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Parameter 

Earth 

Mars - current 

Mars-paleo 

Surface Pressure, P(0) 

1 .0 atm 

.01 atm 

1.0 atm 

Scale Height, H 

7 km 

8 km 

8 km 

Bulk Gas Constant, 

287 J K' 1 kg 1 

191 J K' 1 kg 1 

191 JK' 1 kg 1 

Bulk Molecular Weight, mw 

28.2 g/mol 

43.49 g/mol 

43.49 g/mol 

Gravity, g 

9.8 m s' 2 

3.7 m s' 2 

3.7 ms' 2 

Specific Heat Ratio, y = Cp/C v 

1.4 

1.3 

1.3 


Table 3. Boundary conditions used in examples. 


Parameter 

Earth 

Mars-current 

Mars-paleo 

Radius at Source (rj 

200 m 

1000 m 

100 m 

Velocity at Source (u„) 

300 m/s 

230 m/s 

230 m/s 

Temperature at Source(0 o ) 

1000 K 

1000 K 

1000 K 

Gas Mass Fraction of Magma(n 0 ) 

5% 

1% 

5% 

Mass Eruption Rate (p 0 u 0 rrr 0 2 ) 

1 .6 x 10 8 kg/s 

1.6 x 10 8 kg/s 

1.6 x 10 8 kg/s 

Maximum Predicted Height 

21 km 

133 km 

31 km 
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Figure Captions 

Figure 1 Atmospheric temperatures as a function of altitude assumed in model calculations. 

Temperatures for Earth are for the US Standard Atmosphere. Temperatures for Mars 
and Paleo-Mars are a simplification of those reported by Sief and Kirk [1977] 

Figure 2. Atmospheric pressures as a function of altitude assumed in model calculations. 

Pressures for Earth are for the US Standard Atmosphere. Pressures for Mars and 
Paleo-Mars assumed surface pressures of 1 .0 x 10 3 Pa and 1 .0 x 10 s Pa, respectively. 
Both profiles assumed a scale height of 8 km. 

Figure 3 Mean free path as a function of altitude for Earth, Mars and Paleo-Mars. Also noted 
are the maximum altitude of the Pinatubo eruption column, the maximum altitude at 
which the U2 aircraft can fly, and the maximum observed altitudes for atomic bomb 
tests and firestorms. 

Figure 4. Local speed of sound as a function of altitude for Earth, Mars and Paleo-Mars. Note 
that current and early Mars have similar profiles for speed of sound resulting from the 
dependence on temperature (we have assumed the same temperature profile for both 
cases). Also note that the curves stop at some maximum altitude, above which the 
speed of sound is not defined. 

Figure 5. Bulk vertical velocity as a function of altitude for the Mars boundary conditions given 
in Table 3. The local speed of sound is also shown for comparison. 

Figure 6. Solution for plume radius as a function of altitude for the Mars boundary conditions 
given in Table 3. Note that the radii shown are in kilometers. The full solution 
corresponding to the velocity solution in Figure 5 is shown. 

Figure 7. Solution for plume radius as a function of altitude for the Mars boundary conditions 
given in Table 3. This is the same solution as that shown in Figure 6, but focused only 
on radii between 0 and 120 km. 

Figure 8. Comparison of radial expansion and vertical velocities for the Mars boundary 
conditions given in Table 3. Local speed of sound is also shown. Note that the 
assumption of small radial velocity relative to vertical is violated before the plume 
reaches 1 0 km altitude. 

Figure 9. Comparison of ambient atmospheric density and bulk plume density, as functions of 
altitude, for the Mars boundary conditions given in Table 3. 
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PAHOEHOE TRANSPORT AS A CORRELATED 
RANDOM WALK 

S. M. Baloga and L. S. Glaze, Proxemy Research, 20528 Farcroft Ln, Laytonsville, 
MD 20882, 301-313-0026, steve@proxemy.com 


Submitted to JGR7Solid Earth, December 28, 2001 

Abstract. Transport and emplacement of pahoehoe lava is often dominated by random influences, 
such as budding, small-scale variations in topography, and variations in how the skin cools, stretches, 
and cracks. We develop a random-walk model to describe the transport of small pahoehoe lobes 
when random effects are a significant influence on emplacement. We derive a stochastic formalism 
for random movements on a two-dimensional grid. Our formulation features correlation and bias that 
determine the degree of randomness in the random walk of each parcel of lava. We compare 
qualitative and quantitative predictions of the random walk with observations and data for pahoehoe 
lobes in Hawaii. It is found that a correlated random walk quantitatively describes lobe thickness 
transects whereas the uncorrelated random walk only reproduces the meandering of small 
monofilaments of a few toes. Correlation expresses momentum effects and a chain of memory from 
the front to the source of lava supply. We show that the correlated random walk features a channel- 
forming regime that decays with distance from the source. The theory can also be interpreted in the 
context of a time-dependent planar flow. When many parcels of lava participate, the overall shape 
of a one-dimensional lobe has the appearance of a continuum wave. The types of processes contained 
in the random walk model explain the multitude of concurrent or overriding features found in broad 
fields of toes. The basic formalism can thus be extended to simulate large fields of emplaced 
pahoehoe in either terrestrial or planetary settings. 


1. Introduction 

Broad pahoehoe flow fields, such as 
shown in Figure 1, are generally complex 
collections of much smaller, sometimes 
amorphous features. Such pahoehoe flow 
fields are often built up by the emplacement of 
these small features at many different 
locations sequentially or simultaneously. 
They may be interspersed with a’a lobes in a 
complex time-sequencing of emplacements 
[e.g., Nichols, 1936; Walker, 1972; Swanson, 
1973; Holcomb, 1976], The fact that flow 
fields are not comprised of a handful of 
simple discrete flow units, like well-defined 
and isolated a'a lobes, makes the emplacement 
process extremely difficult to model. Thus, to 


better understand the flow fields as a whole, 
we begin by investigating the emplacement of 
typical small features observed within flow 
fields. The ultimate objective of this work is 
to eventually simulate large flow fields for 
comparison with image and topography data 
for lava flows on Earth and other terrestrial 
planets. 

When one observes the emplacement 
of pahoehoe toes, amorphous fields of toes, 
and certain kinds of lobes, it becomes 
immediately evident that a host of random 
processes influence emplacement. These 
random influences may include topographic 
variations at a variety of scales and small- 
scale random processes associated with the 
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cooling, cracking, and stretching of the outer 
skin. There are also local random effects 
associated with budding, channel, and lobe 
formation, interactions with active or 
emplaced parcels, and fluctuations in the 
supply and distribution of lava. 

There are, of course, systematic and 
deterministic influences on pahoehoe 
emplacement as well. While fluid, a parcel of 
lava will tend to flow downslope. While in 
motion, a fluid element carries momentum 
and is subject to forces from neighboring 
elements or boundaries. One would hope to 
model pahoehoe emplacement with a 
generalized Navier-Stokes equation, 
supplemented by governing equations for 
other processes such as heat transfer, phase 
changes, etc. However, we cannot apply such 
deterministic equations without modeling all 
of the small scale transport processes, 
topographic influences, and local fluctuations 
in lava supply and boundary conditions 
separately. In principal this may be possible, 
but it would be impractical to attempt a 
completely deterministic description of the 
complexity observed in the field. A stochastic 
formalism may be employed to describe both 
the deterministic and random influences that 
act simultaneously on the lava parcels. 

We present here such a stochastic 
formalism that treats the transport of a parcel 
of lava as a random walk. But, what exactly 
does the term “random walk” mean? The 
classical approach to the development of a 
random walk treats a particle beginning at 
some initial location, then moving to all 
subsequent locations based on some 
prescribed set of probability laws. At any step, 
the movement of the particle to the next 
position is a random variable. On any 
particular step, or for any sequence of steps, 
the outcome cannot be predicted with 


complete certainty. All the movements of the 
particle have an element of chance that is 
described by probability concepts, rather than 
exact, deterministic ones. 

There are many different types of 
random walks (e.g., classical uncorrelated 
Brownian motion, correlated, biased, Levy- 
flights, and fractal generalizations [see, for 
example: Chandrasekhar, 1943; Wax, 1954; 
Zauderer, 1983; Klafter et al., 1996; and 
references cited therein]). Some types of 
random walks are trivial to prescribe and 
apply to real physical processes. There are 
also cases where even the mathematical 
formalism defied exposition for almost fifty 
years [Goldstein, 1951; Zauderer, 1983]. 
Even though random effects tend to confuse 
initial observations in the field, it is clear that 
there are different styles of pahoehoe 
emplacement. At one extreme, a simple linear 
filament composed of a few toes may grow at 
the front, extending the filament by one unit. 
Then the process may be repeated. Figure 2 
shows such an example. At the other extreme 
are the highly complex fields that may contain 
toes, lobes of toes, small sheets and channels, 
and other forms that are even difficult to 
name. In complex fields, toes, lobes, sheets, 
and channels can form concurrently and may 
override earlier features [Nichols, 1936; 
Walker, 1972; Swanson, 1973; Holcomb, 
1976; Jurado-Chichay and Rowland, 1995; 
Crown and Baloga, 1999], 

These different regimes of pahoehoe 
emplacement suggest that pahoehoe 
emplacement could be associated with 
different types of random walks. We 
formulate a very general random walk model 
for the movement of parcels of lava on a fixed 
spatial grid. This model determines the 
probabilistic distribution on the grid for 
different types of random walks. However, 
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one cannot determine the type of random walk 
applicable to different emplacement regimes 
from a priori theoretical considerations. This 
must be done by field observations and 
data/theory comparisons. Our model does not 
address the small-scale physics of budding 
and toe formation, but rather, focuses on the 
large scale movements of an ensemble of 
parcels. 

The main step in any random walk 
formulation is to determine the probability 
distribution for the locations of a parcel or 
particle in space or time. Once this is done for 
a given random walk, all that can be known 
about single or multiple particle movements 
(e.g., means, standard deviations, moments) 
can be determined by computation. When 
multiple particles are released from the origin, 
this probability distribution determines the 
concentration of particles at any location in 
space or time. In our application, this 
probability distribution represents the 
thickness of a pahoehoe flow lobe. 

Before presenting the random walk 
model, we first describe some specific 
pahoehoe features found at typical locations in 
Hawaii. We have studied features where the 
influence of random effects was at least as 
important as the deterministic ones. Given the 
wide variety of emplacement behaviors, we 
carefully define what we expect the random 
walk model to accomplish in the subsequent 
section. Next we formulate the emplacement 
of pahoehoe as a random walk characterized 
by two parameters. One parameter determines 
the degree of bias in the random walk, the 
other is a measure of the correlation between 
steps. Readers may be familiar with the 
simplest form of a random walk, uncorrelated 
Brownian motion [e.g., Chandrasekhar, 1943]. 
Our derivation is designed to contain 
uncorrelated Brownian motion as one limiting 
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form. Our mathematical formulation extends 
this random walk so that a bias could be 
included that might serve to represent the 
tendency for movement in the downslope 
direction More importantly, our extension 
also includes the difficult concept of 
correlation among the steps of a random walk. 
This reflects a chain of memory in the 
movements of a parcel due to momentum and 
the influences of neighboring parcels. A key 
aspect of our general formulation is that it is 
set up recursively so that it can be applied 
repetitively to any number of steps or fluid 
elements. 

For comparison, we have made 
observations and collected dimensional data 
on small pahoehoe lobes from many sites in 
Hawaii. We compare the qualitative 
predictions of the model for different types of 
random walks to what is seen in the field. We 
also compare measured cross-sections with 
theoretical predictions. It is shown that the 
correlated random walk describes a host of 
important features found in the field that 
cannot be explained by an uncorrelated 
random walk. As noted in the discussion, the 
theory can also be interpreted in the context of 
a time-dependent two-dimensional planar 
flow. When many parcels of lava participate, 
the overall time-dependent shape of a two- 
dimensional lobes has the appearance of a 
continuum wave. In conclusion, we identify 
some possibilities for future detailed 
data/theory comparisons and extensions to 
simulations of large pahoehoe fields. 

2. Field Observations 

We begin by summarizing our 
observations of active and inactive pahoehoe 
lava flows in Hawaii. The intent of our initial 
observations was to identify the qualitative 
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character of random effects on pahoehoe 
emplacement. These observations were 
expected to guide the development of a model 
of emplacement that would highlight the role 
of random effects. We have observed active 
toe emplacement at the current eruption of 
Pu’u O’o at Kilauea and numerous small 
inactive features on older flow units, including 
the 1969-1974 eruption of Mauna Ulu at 
Kilauea [Swanson, 1973; Holcomb, 1976], 
and prehistoric flows from Mauna Loa and 
Hualalai. We have looked at a range of small 
mappable pahoehoe units, comprised of 
apparently gas-rich to fully-degassed lava in a 
variety of settings, including steep and 
shallow slopes and different substrates 
[Wilmoth and Walker, 1993; Crown and 
Baloga, 1999]. 

For the purpose of developing a 
stochastic model, we have focused on a few 
simple features that permit us to describe the 
nature of the random effect. These features 
include isolated filaments of toes, monolayers 
of toes and similar mappable units [Crown 
and Baloga, 1999], and small lobes that were 
isolated, in pre-existing channels, or in contact 
with other lobes. Simple sheets of pahoehoe 
that appeared to have flowed as a single 
continuum unit according to classical viscous 
dynamics were not of interest for this study. 
Neither do we specifically address other 
formations, such as squeeze-ups from rifted 
tumuli, or overflows from channels [e.g., 
Jurado-Chichay and Rowland, 1995]. 

To develop precise objectives for a 
stochastic model, we must be clear about the 
feature definitions under study and their 
characteristics relative to the model. Our 
descriptions of units and observations are as 
follows: 

Isolated filaments and monolayers: 
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Isolated filaments of a small number of 
toes, typically a few to about a dozen, are 
common in Hawaii (Figure 2). Often such a 
filament advances by budding from the 
forward toe, but toes may also form toward 
the side. There is random meandering as the 
filament advances but there is usually an 
overall directionality to the orientations of the 
toes. When unconfined by pre-existing 
topography, the thickness of such deposits is 
generally equivalent to one or two toes. It is 
not always possible to tell where the next toe 
will form or its orientation, or what the final 
shape of the filament will be. Such filaments 
are often found on steep local slopes (>15 
degrees) and pali. The directionality of the 
filament generally aligns with the downslope 
direction, as would be expected, and the toes 
tend to be distended. However, isolated 
filaments are also found on small virtually flat 
slopes, yet they often retain a directionality. 
This may be due to the initial direction of lava 
from the source or the topographic slope at the 
scale of the toes or filament. In the field, such 
small-scale topography is often difficult or 
impossible to discern. On flat areas, 
emplacement may also take the form a broad 
sheet-like analog of a linear filament. This has 
been called a monolayer [Crown and Baloga, 
1999], Instead of a single strand, there is a 
broad sheet of toes that is typically one or two 
toes thick with a surface of mappable units. 

Isolated and Laterally Confined Lobes: 

We consider a lobe to be isolated if the 
margins were emplaced without obvious 
interference from pre-existing topographic 
influences (Figure 3). Isolated lobes are often 
composed of a few dozen to a few hundred 
parcels of lava. These parcels form toes, buds, 
and a variety of small mappable units when 
they are located at the outer surface of the lobe 
[Crown and Baloga, 1999] Lobes can take the 
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form piles that have transverse profiles with a 
central ridge, a crude bilateral symmetry, and 
a decay in the thickness of the deposit toward 
the margin. The lobate features we discuss 
here typically have a width of several meters 
to as much as ten meters and a length of about 
five meters to as much as a hundred. The outer 
surface of the deposit is composed of toes and 
there are often interior parcels of lava that 
never had the opportunity to form toes. 
Transverse thickness profiles initially seemed 
to be reminiscent of bell-shaped Gaussian 
cross-sections, modified by the the discrete 
nature of the toes themselves and the natural 
variability of the toe formation process. At a 
small scale, the lobe emplacement process is 
dominated by the formation of individual toes. 
However, at the scale of the lobe itself, the 
advancing motion may often resemble a 
continuum flow. 

Where lobes are laterally confined, the 
cross-section shape tends to have a relatively 
flat centroid and a relatively rapid decay 
toward the margins. As with isolated lobes, 
this is sometimes reminiscent of a bell-shaped 
curve, although the toe formation process is 
not allowed to continue past the margins. 
Figure 4 shows a confined lobe where the 
lateral extent of the lobe was influenced by 
pre-existing topography. Unfortunately, 
“lateral confinement” is not always easy to 
distinguish in the field, as we initially 
preconceived. Based on a general qualitative 
survey of dozens of such features, it appears 
that the flow itself often provides its own 
lateral confinement through the emplacement 
of toes. This suggests that lateral self- 
confinement is as important, or perhaps even 
more important than the ambient topography 
in determining the cross-sectional shape of the 
deposit. Although more data is needed, we do 
not believe there is a systematic difference 
between the cross-sectional shapes of lobes 


that are confined and those that are not. 

Broad amorphous fields: 

Broad, complex flow fields result 
when multiple modes of lava emplacement 
occur within the same area or are 
superimposed on each other. Such fields may 
be composed of hundreds or thousands of 
toes, isolated and confined lobes, sheets, 
channels, and other features. When a broad 
field is still active, one can observe a toe form 
in one location while other toes intermittently 
form tens of meters away. As one area 
becomes inactive, others may become active 
at seemingly random locations in the field. 
This suggests very complicated spatially 
random lava transport patterns within the flow 
field (Figure 1). 

3. Objectives for the stochastic model 

If we must abandon a detailed 
deterministic approach to the emplacement 
process, one might reasonably ask, “What can 
we expect to accomplish with a stochastic 
model?” Based on the observations above, we 
have selected a limited set of tractable, 
practical objectives. The random walk model 
is expected to: 

* Provide at least a qualitative 
description of the emplacement of 
monofilaments of small numbers of 
toes, 

* Describe the spatial spreading rate of 
the margins when there are enough 
parcels to form a thin layer or lobe, 

* Approximate quantitatively the shape 
of transects across lobes, 
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* Retain contact, under some set of 
parameters, with conventional 
descriptions of pahoehoe emplacement 
as a fluid dynamic and heat transfer 
process, 

* Qualitatively describe transitions from 
emplacement regimes that are 
dominated by random motions to 
those that are coherent, organized, and 
subject to classical fluid mechanics, 
and, 

* Suggest what should be measured in 
the field at length-scales comparable 
to toes and at length-scales 
comparable to fl ow lobes. 

From a more altruistic perspective, we 
would like the random walk model to be based 
on sound probability theory. It should produce 
dimensions or features that can be compared 
with field data and contain the possibility of 
relating the parameters of the model to 
conventional physics of mass, momentum, or 
heat transfer. We can of course simulate any 
sequence of random events on a computer 
provided enough decision criteria are built 
into the simulation. We would like the model 
to provide a general means for discriminating 
whether a particular simulation is 
representative or not. 

4. The Correlated Random Walk 
Model 

To formulate and interpret a random 
walk model of lava flow emplacement, we 
must be specific about the definition of a 
“parcel” of lava. A parcel is considered to be 
a fixed volume of lava that is sufficient to 
form a typical toe or some other form of 
discrete mappable unit. While in motion, 


possibly within a much larger ensemble of 
such parcels, it may be a fluid element that 
actually bears no resemblance to its final 
form. We do not suggest that all parcels 
become toes, but merely use such a concept as 
a convenience. In our work here, we treat 
parcels as having a single volume. Statistical 
variations in parcel dimensions in Hawaii and 
Idaho are given in Crown and Baloga [1999], 
but such differences are not treated here. 

Based on field observations, it is clear 
that the random walk formalism must include 
the possibility of a bias. Such a bias 
represents a tendency for motion to occur in a 
preferred direction. For example, an isolated 
parcel may tend move in the direction of the 
downslope gradient. Also, parcels may carry 
momentum and be influenced by prior 
motions and source conditions. These 
processes are more closely associated with the 
concept of a correlation between movements 
rather than a bias. So, the model should also 
incorporate the possibilities of perfect 
correlation, its absence, and the full range of 
intermediate cases. To be useful for 
simulating large flow fields, or planetary scale 
settings, the theory must be applicable to a 
full range of scales from a few toes to 
thousands of parcels. We desire a general 
analytic formalism that would allow us to 
relate the model to conventional fluid 
dynamics. Finally, for practical reasons, the 
formulation must avoid having to manually 
derive the probability distribution for each 
step for different degrees of correlation or 
bias. 

Our formulation of the correlated 
random walk uses two parameters to describe 
different types of random walks on the grid 
shown in Figure 5. See Table 1 for a 
definition of mathematical symbols. The first 
parameter a is the probability of taking an 
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upward (+y) step from any location when the 
choice is random. Thus, l-a is the probability 
of taking a step in the downward (-y) 
direction. When a = 0.5, the random walk is 
said to be unbiased. 

We also use a parameter p that 
provides a measure of the correlation in the 
random walk and the amount of randomness. 
When p =0, the random walk is said to be 
uncorrelated. This means that the tendency of 
move upward or downward at any step is 
governed solely by a and is completely 
independent of any choices made on prior 
steps. When p =1, there is no randomness 
except for the first step. We refer to this case 
as perfect correlation. Once a parcel starts in 
a direction, it continues identically in that 
same direction. In general, for arbitrary 
combinations of the parameters, a random 
walker can be at any y- value from -N to +N 
after N steps. We seek a general formulation 
for any combination of a and p values. 

The objective is now to determine the 
explicit form of P(y, N), the probability that a 
parcel is located at a position y on the grid 
after N steps have been taken. We use 

a recursive matrix formulation that determines 
the probability distribution after N steps based 
on the distribution after N-l steps: 

'0, AO] = [A (N)][P(y, AM)], N >(1) 

The key to such a formulation is finding the 
form of the [N-l x N] matrix [A (N)]. Once 
[A (N)] is determined explicitly, eq 1 can be 
used to obtain the probabilities at any step for 
any mixture of bias and correlation. The 
recursive process is started by specifying [P(y, 
N= 1 )] with the two elements P(y=+1, N=l) 
=a and P(y=-1, N=l) = l-a. All subsequent 
probabilities are then found by repeated 
applications of eq 1. 


Determining the explicit form of [A 
(N)] was very challenging because the 
correlations originate directly from two prior 
steps, and indirectly from all prior steps, while 
the recursion relation is based only on the 
immediately prior probabilities. As part of 
this determination, we must also ensure that 
(1) the probabilities retain normalization: 
£P(y;N)=l at each step, (2) for p=0, we 
recover the classical binomial distribution 
describing an uncorrelated random walk with 
an arbitrary bias, and (3) when p=l the initial 
motions are retained identically. 

For P(y, N) to be a rigorous 
probability, it must satisfy the normalization 
requirement (1). This imposes constraints on 
the elements of the matrix [A], If we require 

E A, , (A) = 1 ^2) 

J 

i.e., the elements of each column sum to 1, 
then P(y, N) will be normalized if P(y, N-l) is 
normalized. This is proven by summing all 
the P(y, N) on the left side of eq 1 , writing out 
the right hand side explicitly, and using eq 2. 
Thus, 

E Piy, AO = 

iV ■. y \ utM o,„ ),A 

E P(y,N-\) 

.V « lilt eirn < otid ) i N I 

So, once we assure that the probabilities of the 
first step are normalized, the probabilities of 
all subsequent steps will be guaranteed if the 
elements of [A] satisfy eq 2. 

The requirement (2) to reproduce an 
uncorrelated random walk means the elements 
of [A] must make: 
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P(y; A0= 

W{ k,N,a) = W((N+y)/2; N;a) = 

M , „ „ , ( 4 ) 

a (1-a) ‘ 

[(A + y)/2]![(A-.y)/2]! 

when p = 0. Here W is the binomial 
distribution for k= (y+N)/2 successes out of N 
Bernoulli trials with success probability a. It 
is interesting to recall here that the limit of the 
binomial distribution as N becomes large is 
the classical Gaussian or normal distribution 
that describes bell-shaped curves. Our 
requirement (3) to retain initial motions 
identically when p = 1 means that [A] must 
make: P(y=+N, N)=a and P(y=-N, N)= 1-a for 
the case of perfect correlation. 

We have derived the general form of 
[A (N)] by analyzing the grid in Figure 5 and 
the three requirements above. Now we 
indicate, by example, how the matrix elements 
were obtained. Correlation is, in effect, 
memory of the choice made by a parcel in the 
prior step to travel upward or downward. 
Thus, it first becomes evident on step N = 2. 
So, the first task is to obtain A (N=2) from the 
three requirements and the recursion relation: 

[P(y, 2)] = [A(2)][P(y, 1)] (5) 

The derivation of the six elements of 
[A (2)] is easily illustrated by the tabular form: 



-12 

?(:,!; = it 


4 2 

?(-i ,i) = !■■& 


•ft* 



If we examine Figure 5, we see that the parcel 
cannot travel to y=+2 if the first step was 
downward or y = -2 if it originally went 
upward. Consequently A, 2 = A 3I = 0. If a 


parcel travels to y=+2, it must have 
probability a 2 when p = 0 and probability a 
when p = 1. This can be obtained by setting 
A, | — p+a(l-p) = p ++ . Here we introduce the 
superscripts because these elements will 
appear elsewhere in the general matrix for 
[A(N)]. To satisfy the normalization 
condition of requirement (1), this further 
implies A 21 = 1 - p - a(l-p) s p + * . By similar 
analysis of y = -2, we must set A 32 = 1 - a ( 1 - 
p) = p" to obtain the probability (1-a) 2 when 
p = 0 and probability 1-a when p = 1 . Then 
normalization requires A 22 = a(l-p) = p' + . 

All six of the elements of [A (N=2)] 
have thus been determined explicitly in terms 
of the two parameters a and p. It is easy to 
verify that the P(y, 2) sum to 1 for any 
combination of a and p, the binomial 
distribution is obtained for p =0, and 
P(y=+2, 2)=a and P(y=-2, 2) = 1-a when p = 
1. 

The same seemingly simple process 
can be used to obtain [A (N=3)] and [A 
(N=4)]. However, as the reader may wish to 
verify, the derivation rapidly becomes much 
more complex. Unforeseen problems arise 
with the edge effects, counting the central 
paths correctly, and retaining normalization of 
P(y=3, 3) and P(y=4, 4). Only when N > 4, do 
we appreciate the full difficulty of the 
problem. How this occurs can be seen from 
Figure 5. If we take an arbitrary location on 
the grid, we see that there are, in general, four 
possible locations that can result in the two 
subsequent steps. Consequently, if we use the 
mirror image of such a wedge, we impinge on 
the boundary in different ways coming 
backwards from step 4 to step 2. 

We have found, nevertheless, our main 
theoretical result to be the general explicit 
form of [A (N)] as: 
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p ++ 

0 

0 

- 


0 

P T ‘ 

p- T++ 

0 




0 

p- + - 

p- +++ 




0 

0 

p-- + - 




1 










p- + - 

p- + 

0 




0 

p” 


where the probabilities p' +++ and p" + ‘ are 
defined in Table 1. We note that this matrix 
always has dimensions of [N-l x N], and that 
only the diagonal elements shown contain 
non-zero values. 

By using this result in eq 1, we can 
readily generate the probability that a parcel 
will be located at any admissible y-value after 
N steps. It is straightforward to show that all 
three requirements are met in the general case 
by using the explicit elements for the [A (N)] 
from Table 1 . This general form for [A (N)] is 
the main result of the theoretical formulation. 
With it, we can obtain the distribution of an 
ensemble of parcels for any desired 
combination of bias and correlation. When 
many random walkers are considered, the 
thickness of a lobe is considered to be 
proportional to the probability that a single 
random walker is at location y at step N. 

5. Applications 

Monofilaments: 

The simplest application of the random 
walk formalism is the advance of a 
monofilament. We show this trivial case to 
illustrate future possibilities for more detailed 
simulations. Here our goal is to model the 


meandering of the strand. We do not consider 
here the possibility of budding at locations 
other than the toe, although this clearly occurs 
in the field. As a toe forms at the front of its 
predecessor, it is considered to take a single 
step in the x-direction. Simultaneously the 
parcel is also considered to take a random step 
in either the plus or minus y- direction. A 
filament of N toes is thus formed by a 
sequence of N random steps. The process of 
simulating a monofilament corresponds to 
generating a single path in the network shown 
in Figure 5. 

As an example, we suppose that an x, 
y coordinate system is laid out in an arbitrary 
direction and that 8 toes are to be emplaced 
along the path of the random walk. We 
consider the random walk to be uncorrelated 
but biased with a= 0.75. This means the 
probability of taking a step in the +y direction 
is 0.75 and the probability of taking a step in 
the -y direction is 0.25 for each step of the 
sequence. On the average the monofilament 
has a constant drift along our coordinate 
system as long as a is constant. In our 
example, the average step is then <y> = 0.5N, 
and thus, after 8 steps we expect the tip of the 
filament to be at y=+4. 

We have simulated such paths with a 
random number generator by setting 

Ay = e V i ( 6 ) 

where the e are taken from a Bernoulli 
distribution with <y>=0.5. It should be note 
that many random number generators produce 
the Bernoulli distribution P(z) for z £ [0, 1]. 
To convert this to our y coordinates, we use 
the transformation 
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Results of these simulations are shown 
in Figure 6. The straight line shows the 
average path. Other paths correspond to 
different simulations with the same 
parameters. It is very clear that any particular 
sequence may have a very different shape 
from the average. The randomness in the toe 
orientations is evident, as well as the influence 
of the bias in each step. Physically the 
randomness in the orientation could be related 
to the location of the budding site, while the 
bias could be due to the overall downslope 
direction with respect to our coordinate 
system. 

This trivial example reflects 
qualitatively the monofilament meandering 
found in the field. Field studies of such 
orientations could be performed. However, we 
anticipate that little would be gained because 
this style of emplacement is dominated by the 
detailed microphysics of the budding process 
and the stochastic formulation is directed 
towards larger ensembles. We leave it to the 
interested reader to determine the variance or 
standard deviation of this application and its 
conversion into the distribution of orientation 
angles. The interested reader may also wish to 
tackle similar simulations with the correlated 
random walk 

At, = 0 -p) e +pAy ( V / (8) 

and investigate the influence of correlation on 
the angular distributions of orientations. 
Clearly, the discrete steps of our trivial 
example can be extended to continuous 
distributions as given in Crown and Baloga 
[1999], 


Lobes: 

An isolated lobe is produced by the 
model when there are enough parcels to form 
a distribution P(y, N) at all accessible y 
values after N steps. Thus the lobe will have 
relief and some shape in cross-section. When 
correlation is absent, the result is the cross- 
section given by eq 4. Figure 7 shows typical 
theoretical transects for the uncorrelated 
random walk at three different advancing 
steps of lobe emplacement. 

We have visited numerous sites at 
Kilauea, Mauna Ulu, Haulalai and various 
other locations in Hawaii, where lobes are 
readily found on the flow fields. We found 
that isolated lobes sometimes have a mild 
medial peak with an approximately symmetric 
decay toward the boundaries. Sometimes 
transects are reminiscent of the bell-shaped 
curve. Our initial qualitative inspections of 
isolated lobes suggested that the emplacement 
process might simply be an uncorrelated 
random walk. These qualitative observations 
also suggested that a more quantitative data- 
theory comparison was desirable. 

Figures 8a and b show transverse 
thickness profiles for two typical lobes at 
Hualalai. The measurements were obtained 
by stretching a leveled rope above the lobe, 
perpendicular to the general direction of flow, 
and measuring the relief to the nearest 
centimeter at evenly spaced intervals. Because 
the surfaces of the lobes are irregular, a 
replicate transect was measured at the same 
station in each case to indicate the natural 
variability. 

We have collected about two dozen 
such transects to see if the hypothesis of an 
uncorrelated random walk was consistent with 
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the data (Figure 8). These figures do indeed 
show a crude connection with the measured 
transects and the theoretical cross-section 
associated with an uncorrelated random walk. 
However, the most striking difference one 
observes in the field is that the lobes are 
generally flatter on the top than one might 
attribute to a Gaussian. Although the random 
walk approach seemed promising, the most 
elementary uncorrelated form does not 
provide a good quantitative or qualitative 
comparison with field transects. 

It is not completely surprising that the 
applicability of the uncorrelated random walk 
to pahoehoe transport is limited. First, this 
classical random walk minimizes the 
interaction between parcels during transport. 
Thus, the limited influence of neighboring 
parcels precludes the application of 
conventional fluid dynamics to the entire 
active lobe. If the uncorrelated walk was 
always valid, we would have to abandon all 
conventional viscous fluid dynamic models 
except for a parcel-by-parcel applications. A 
second problem with the uncorrelated random 
walk became evident from computer 
simulations of many uncorrelated random 
walks. These simulations superimposed many 
monofilaments. They reproduced broad, 
amorphous flow fields of many thousands of 
parcels, but were far too diffuse to represent 
solitary or multiple lobes. 

We have addressed these issues by 
considering the correlated random walk. With 
a correlated random walk, some “memory” of 
what happened to a parcel at prior steps is 
retained. Correlation is a manifestation of the 
momentum of the flow, inertial effects, and 
fluid dynamic pressure when there are many 
contiguous parcels in a flow lobe. 

Figures 8a and b also show the 


theoretical cross-sections when correlation is 
present. The flatness of the shape that 
corresponds to what is typically observed in 
the field is obtained when the correlation 
parameter p is in the vicinity of 0.25 to 0.30. 
For lesser values, the cross-sectional shape is 
more leptokurtic. Although the flatter cross- 
sections are typical, some lobes do show a 
modest medial peakedness. 

Figure 9 shows the distribution of 
pahoehoe parcels for a higher level of the 
correlation parameter p. The profiles are 
shown at subsequent time steps from 5 to 20 
with p =0.55. Two noteworthy features are 
evident. First, correlation at this level 
produces noticeable lateral margins in the 
upstream region. The shape of the margins 
clearly indicates a tendency of the flow to 
channelize as a direct result of the correlation. 
Such an effect has been documented [Crown 
and Baloga, 1999] in field studies of the 
transitional regime from discrete pahoehoe 
toes to larger units. Such flow margins cannot 
be produced by an uncorrelated random walk 
without lateral topographic confinement. 
Second, the influence of the correlation 
eventually decays for sufficiently long 
distances from the source. After 20 steps, all 
influence of correlation has been lost and the 
transect cannot be distinguished from a 
normal distribution. Qualitatively, this is what 
is often observed in the field. At the front of 
a small lobe, there can be budding and the 
random formation of new units while lava 
transport upstream is more coherent. 

With an uncorrelated random walk, all 
the influences on the formation of a toe or the 
transport of a parcel are local. The 
uncorrelated random walk describes a number 
of qualitative features observed in the field, 
such as the meandering of isolated filaments 
and the budding of toes. This model seems to 
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be adequate for a small number of toes, say 
less than about a dozen, that occupy a single 
path on the grid. The uncorrelated random 
walk predicts a medial thickening of the lobe 
in cross-section with a gradual decay to the 
margins, which is not the typical case 
observed in the field. 

The influence of correlation becomes 
a significant effect in the transport of 
pahoehoe as the number of parcels increases 
and large-scale fluid dynamic effects compete 
with the random processes. As with the 
classical correlated random walk [Zauderer, 
1983], our formalism also shows that any 
correlations due to pressure, momentum, and 
inertia must eventually decay. In the 
continuum limit for very large numbers of 
parcels, correlated random walks engender 
wave equations [Zauderer, 1983; Goldstein, 
1951], These wave equations, in turn, can 
naturally describe surface wave structures and 
lobate features at the fronts and margins, as is 
often seen in terrestrial and planetary lava 
flows. From classical statistical physics [cf 
Zauderer, 1983 and refs therein], correlation 
can be associated directly with the terms used 
for viscous fluid flow. It also provides a 
“chain of memory” between downstream 
parcels and the source of lava supply. In 
effect, the correlated random walk provides a 
bridge between random effects at the toe or 
parcel-scale and large-scale continuum fluid 
dynamics. 

6. Discussion 

It must be admitted that a possibility 
for obtaining better agreement of the 
uncorrelated random walk with field 
observations is by lateral confinement at the 
margins. Whether by pre-existing topography, 
or perhaps by the cooling of the outermost 


lateral boundaries, the effect would be to 
thicken the margins, thus flattening the cross- 
section. For completeness, we now determine 
the transect shapes that are produced by a 
laterally confined uncorrelated random walk. 

We suppose there is some confining 
boundary on one side at y 0 . When 
uncorrelated random walkers encounter this 
boundary, they are reflected back into the 
admissible parts of the network and resume 
their random walk. Chandrasekar [1943], 
with reference to Smoluchowski [1915], uses 
a clever analysis of the image of paths beyond 
the reflecting boundary to show that, for a 
single reflecting plane at +y Q , 

P(y\N)= 

W((N+y)/ 2; N;a) + P(2_y -y; N) ^ 
\{N+y)!2] N\ a) + W((N+2y\ -y)f 2; A 


It is straightforward both conceptually and 
analytically to verify that the probabilities 
given by eq 9 retain normalization. 

We have used this approach with two 
reflecting boundaries at +y 0 and -y 0 to 
represent the lateral confinement of lobes. To 
use this formula, one must be careful about a 
number of details. For example, on the array 
we have specified, the sites that are accessible 
depends on whether N is even or odd. For 
N=even, y must be even, and the reflecting 
boundary must be from the prior step. This 
means the boundaries for lateral confinement 
must be at odd locations. For N=odd, these 
considerations reverse. The lateral 
confinement of the random walk has the effect 
of trimming off the tails by shifting the parcels 
to admissible interior parts of the distribution. 
Due to normalization, the effect is to steepen 
the sides of the distribution, thus flattening the 
overall shape. Figures 9a and b show two 
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transects that more closely resemble what is 
observed in the field. 

At present, we cannot preclude the 
laterally confined uncorrelated random walk 
as a viable model for explaining the cross- 
sectional shape of pahoehoe lobes. However, 
there are three points that argue against it. 
First, the transects we have shown in Figures 
8a and b represent typical cases. For the 
confined uncorrelated random walk to 
reproduce such shapes, approximately the 
same amount must be trimmed from the tails 
in actual lobes. As yet, we know of no 
physical mechanism that can cause this to be 
the usual behavior observed. Second, we 
would expect to find in the field a full 
spectrum of confinements, ranging from the 
Gaussian shape to extreme lateral 
confinements. This is not consistent with out 
observations, even when a pahoehoe lobe is 
found within embanking channels and 
between pre-existing flows. Third, it is 
physically appealing that, when enough 
parcels are present, and the random influences 
are relatively small, the lobe should be 
governing by classical continuum viscous 
dynamics. Such a limit cannot be attained by 
the uncorrelated random walk. 


lobe or the flow field. 

Mathematically, such a limit is 
obtained in the following way. We suppose 
that in a small quantity of time A, a parcel 
takes a step of length 6 along the x-axis 
according to the Bernoulli distribution. We 
consider three bins at x-S, x, and x+S, and ask 
how the probability P(x,t) (or equivalently the 
total number of parcels) at x changes in the 
time A. The change is given by 

P(x, t+ A) - P(x,t) = 
aP(x-5,t))+(\-a)P(x+b,t) - (10) 

(1 a)P(x,t) - a P(x,l) 


If we expand this equation in a Taylor series 
that is first-order in time and second-order in 
space, we obtain, after taking the limits of the 
small quantities. 


— + 

dt dx Qt 1 


( 11 ) 


where 


V 


( l-2a ) — , D 

A 


6 ' 


2 A 


( 12 ) 


We now turn our attention to a time- 
dependent formulation of the random walk 
model for a 2-dimensional flow. Our purpose 
here is to illustrate how the parcel-scale 
considerations can be used to develop macro- 
scale governing partial differential equations. 
For simplicity, we use only the uncorrelated 
random walk. We now employ a continuum 
limit to show what happens when the number 
of parcels becomes infinitely large. In such a 
limit, the small-scale motions of individual 
parcels are considered beneath the resolution 
of the model and are indistinguishable. Thus, 
we observe only the overall evolution of the 


This is the classical diffusion equation 
with a kinematic transport term, where the 
diffusion coefficient D and the “wave 
velocity”, V, are given in terms of the random 
walk parameters A and 
6. We can readily solve this equation for a 
constant flux at x=0. The time-dependent 
solutions are shown in Figure 11. In the 
continuum limit, even for the case of an 
uncorrelated random walk, we reproduce the 
shape of an advancing wave that looks much 
like what we would expect from conventional 
fluid dynamics. It remains for future 
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investigations to develop the correlated 
random walk in such a continuum limit. 

Conceptually, we can now begin to 
interpret or simulate much of what is seen on 
the surface of broad, seemingly amorphous 
pahoehoe flow fields. The model represents 
several regimes of emplacement, including 
small monofilaments and monolayers of toes, 
small lobes that may have channels, and larger 
lobes that advance as a continuum wave. 
Flow fields are often superpositions of these 
different modes of emplacement. We do not 
at present know what controls the formation 
of one style of emplacement versus another. 
The quanitity, location, and timing of the lava 
supply are clearly a major factors. However, 
flow rate, cooling rate, texture, setting and gas 
content may all be important [Rowland and 
Walker, 1987; Rowland and Walker, 1990; 
Keszthelyi and Denlinger; 1996;Wilmoth and 
Walker; 1993; Crown and Baloga, 1999]. 
These, in turn, have there own random 
variations, which remain to be studied 
subsequently. 

7. Conclusions 

The transport and emplacement of 
pahoehoe often features random effects that 
compete with, or dominate, systematic effects. 
When such random effects are important, the 
use of a stochastic approach is necessary. 
This means that we can only predict 
dimensions, shapes, and features in a 
probabilistic sense. Even under identical 
starting conditions, particular realizations may 
differ from case to case. It is impossible to 
predict the details on any particular 
realization. Length and width become random 
variables, as well as any simulation of 
locations, thicknesses, etc. When the 


emplacement can be modeled as a random 
walk, the model provides a means for 
evaluating whether a particular case is 
representative or not. 

The random walk theory we have 
presented embraces different types of random 
walks, biased, uncorrelated, and correlated. 
These three effects can be adjusted in the 
general case by two parameters, a and p. 

The first parameter describes the pre- 
existing slope or other factors that bias the 
orientation of motion. The second parameter 
describes correlation between steps that is 
conceptually associated with momentum, 
inertia, and pressure effects. Small-scale 
pahoehoe features, such as monofilaments of 
toes and small lobes are amenable to 
description as different types of random 
walks. 

We have collected field data intended 
to distinguish the different types of random 
walks for the small pahoehoe features found 
in flow fields. We have shown that an 
uncorrelated classical random walk does not 
completely describe what is typically observed 
in the field, except for the meandering of 
small filaments. Correlation must be added 
to the random walk to describe larger lobe 
cross-sections qualitatively and quantitatively. 
A special case of the uncorrelated random 
walk with lateral confinement can also 
reproduce the rather flat lobe thickness 
transects observed in the field. Although this 
special case cannot be rigorously precluded, a 
host of factors argue against its validity. 

The correlated random walk provides 
an adequate description of lobe thickness 
transects and provides a conceptual bridge 
between random processes and conventional 
viscous fluid dynamics. It can be related to 
momentum transport, hydrodynamical fluid 
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pressure, and inertial effects, all of which 
operate to some extent when many lava 
parcels are contained in the active lobe. The 
natural spreading of a correlated random walk 
clearly has a channel-forming regime that 
decays to uncorrelated behavior far from the 
source, as is often seen in the field. 

We have shown for the uncorrelated 
random walk how to obtain macro-scale 
governing partial differential equations from 
the random model at the lava-parcel scale 
when hundreds or thousands of parcels are 
involved in the emplacement. The time- 
dependent longitudinal transects thus obtained 
resemble the advance of a continuum wave. 
However, the more exciting prospect resides 
with the richness of possibilities, e.g., 
channels and wave-like structures, contained 
in the macro-scale correlated random walk. 

The underlying physical processes and 
factors that determine the type of random walk 
in particular field settings are not yet 
understood. Future studies will focus on this 
connection, with emphasis on the transitions 
for one type of random walk to another and 
the explicit relationships to conventional lava 
transport models. The investigations reported 
here suggest that the random effects 
associated with pahoehoe emplacement are 
indeed amendable to quantitative treatment 
and offer promise for modeling the 
emplacement of planetary flows. 
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Figure Captions 

Figure 1: A typical broad amorphous 

field of toes in Hawaii. This figure 
indicates the predominance of random 
effects on the emplacement process 
and the need for a stochastic approach. 

Figure 2: An isolated filament formed by 

budding of toes. For the most part, the 
filament advances by the budding 
from the forward toe, but toes also 
form toward the side. There is random 
meandering as the filament advances 
but often a general overall 
directionality with randomness in the 
orientations of the toes. 

Figure 3: A typical isolated lobe of 

pahoehoe toes. Note the thickening 
medial tendency and the decay toward 
the laterally unconfined margins. The 
transect shape is reminiscent of a bell 
shaped curve. 

Figure 4: A typical laterally confined 

lobe. Here the lateral extent of the 
lobe was influenced by pre-existing 
topography. 

Figure 5: 2-dimensional grid for the 

random walk. N and y are discrete 
variables representing the downflow 
distance and lateral distance 
respectively. 

Figure 6: Simulations of a random walk 

for a monofilament. The straight line 
shows the average path. Other paths 
correspond to different simulations 
with the same parameters. 
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Figure 7: Typical theoretical transects 

for the uncorrelated random walk at 
three different advancing steps of lobe 
emplacement. 

Figure 8 a, b: Typical examples of thickness 
cross-sections of isolated lobes at 
Hualalai. Each panel shows two 
replicate sets of measurements at the 
same station. Also shown are 
theoretical cross-sections with p-0.28 
as solid line. 

Figure 9: Theoretical cross-sections of lobes 
produced by a correlated random walk 
with p=0.45 at time steps ranging 
from 5 (top curve) to 20. 

Figure 10 : Cross-sections of a laterally 

confined lobe and comparison with 
laterally confined uncorrelated random 
walk. 

Figure 11: The wave formed by the 

advance of uncorrelated random 
walkers in the continuum limit. The 
two waves farthest right at t=75 and 
t= 1 00 have a=0.9. The two waves to 
the left at t=75 and t= 1 00 have 
a=0.75. 
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Table 1. Definition of Mathematical Symbols 


Symbol 

Definition 

y 

Distance from origin/ flow center 

N 

Step number/downflow distance 

0< a< 1 

Single step probability (+y direction) 

0< p< 1 

Correlation parameter 

[A] 

Recursion matrix 

P t+ 

p+a(l-p) 

P + ' 

l-p-a(l-p) 

P + 

a(l-p) 

P" 

1- a(l-p) 

p- + - 

(p- + p + -)/2 

P + " 

(P + +P ++ )/2 

W(k;a, N) 

Binomial distribution for y successes 
in N Bernoulli trials with success 
probability a 
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